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Abstract 

Bead packs of up to 150,000 mono-sized spheres with packing densities ranging from 
0.58 to 0.64 have been studied by means of X-ray Computed Tomography. These 
studies represent the largest and the most accurate description of the structure of 
disordered packings at the grain-scale ever attempted. We investigate the geomet- 
rical structure of such packings looking for signatures of disorder. We discuss ways 
to characterize and classify these systems and the implications that local geometry 
can have on densification dynamics. 
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1 Introduction 



For several centuries, the structures generated by tightly packed mono-sized 
spheres have been studied by scientists as model systems for understanding 
the emergence of order and crystallisation or -vice versa- the appearance of 
disordered and amorphous phases in natural systems [1]. Despite all the ef- 
forts, several questions remain unsolved. One of the most intriguing and chal- 
lenging problems is to understand whether or not the notion of 'ideal' or 
'typical' or 'common' disordered packing configurations has any empirical ba- 
sis [2]. Indeed, disorder does not exclude organization and the positions of 
the packed spheres are locally highly correlated. Experimentally, it is easy to 
observe that when spheres with approximately equal sizes are poured into a 
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Fig. 1. (left) Volume rendering of ~ 150, 000 sphere-pack in a cylindrical container, 
(right) Same image with the topological distances from a given central sphere high- 
lighted in colours (online version). 

jar, squeezed into a bag, shaken or mixed, they form structures with densi- 
ties in a narrow range. No stable packing configurations have been observed 
with densities below 0.56 and no disordered assemblies of spheres have been 
observed with densities above 0.645 [1]. Molecular dynamics simulations show 
that hard sphere systems behave like a gas below p ~ 0.49. Upon compaction, 
liquid-like behaviour is observed up to p ~ 0.55, where crystallisation starts to 
occur. If crystallization is avoided, the system undertakes a 'glass' transition at 
p ~ 0.56 and then it can be compacted up to p ~ 0.645 where no further den- 
sification can be induced [3]. Empirical and simulated evidence suggest that 
something very special might happen in the geometry of the packing at den- 
sities above p ~ 0.56; a process which must terminate below p ~ 0.65. What 
makes the understanding of this process particularly challenging is that there 
are no a priori reasons for the densification process to stop around p ~ 0.64. 
On the contrary, there are plenty of local configurations which are denser than 
this limit and any arrangement of spheres in stacked planar hexagonal closed 
packed layers (the so called Barlow packings [4]) can reach the density of 
p = 7r/V^8 ~ 0.74, as achieved in fee (face-centered cubic) or hep (hexagonal 
closed-packed) crystalline packings. 



2 Experimental Apparatus and Methodology 

In this paper we report the analysis on 6 bead pack samples, all in a cylindrical 
container with an inner diameter of D ~ 55 mm and filled to a height of 
H ~ 75 mm. Specifically: 
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• Two large samples with 150,000 beads and diameters d = 1.00 ± 0.05 mm; 

• Four smaller samples with 35,000 beads and diameters d = 1.59 ±0.05 mm. 

A packing realization is shown in Fig.l. To our knowledge, the present study is 
the first empirical investigation of large packings of monosized spheres. Indeed, 
there are only two previous investigations which use a similar approach [5,6], 
but the first concerns a single test-sample with ~ 2, 000 spheres only, whereas 
the second investigates polydisperse assemblies with ~ 15, 000 spheres. 

Sample Preparation: Low density packings are obtained by placing a stick 
in the middle of the container and pouring the beads into it. Then, the stick 
is removed leaving a packing density p ~ 0.58 [1]. Higher densities, up to p ~ 
0.63, are achieved by gently tapping the sample. The densest sample at p = 
0.64 is obtained by a combined action of gently tapping and compression from 
above (with the upper surface left unconfined at the end of the preparation). 
To reduce boundary effects, the inside of the cylinder have been made rough 
by randomly gluing spheres on the internal surfaces. 

XCT Imaging: A X-ray Computed Tomography apparatus (see Sakellariou 
et al. in this Issue [7]) is used to image the samples. The two large samples 
were analysed by acquiring data sets of 2000 3 voxels with a spatial resolution 
of 0.05 mm; whereas the four small samples were analysed by acquiring data 
sets of 1000 3 voxels with a spatial resolution of 60 pm. After segmentation (see 
Sheppard et al in this issue [8]) the sample data sets are reduced to three- 
dimensional binary images, representing two distinct phases, one associated 
with the spheres and the other with air space. 

Sphere Centres: In order to proceed with the analysis of the geometrical and 
statistical properties the position of all sphere centres are calculated from the 
binary images. Our approach is to find the sphere centres by moving a reference 
sphere (S) throughout the binarised sphere pack (P) and measuring the local 
overlap between S and P. This corresponds to a 3-dimensional convolution: P* 
S. This method is made highly efficient by applying the convolution theorem 
which allow to transform the convolution into a product in Fourier space: 
T\P * S] — jFfPjjFfS'], where T represents the (fast)Fourier Transform. The 
algorithm proceeds in 4 steps: 1) fast Fourier transform of the binary image 
(JF[P]); 2) transform the digitised map of the reference sphere (.F [«!?], chosen 
with a diameter about 10 % smaller than d); 3) perform the direct product 
between these two; 4) inverse-transform of the product: jF~ 1 [J r [P]^ r [S']] = 
P* S. The result is an intensity map of the overlapping between the reference 
sphere and the bead pack, where the voxels closer to the sphere centres have a 
higher intensity. A threshold on the intensity map, locates the groups of voxels 
surrounding the sphere centres. The centre of mass of these grouped voxels is 
a very good estimation of the sphere centres in the pack. 
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0.630* 


0.01 


19,843 
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7.40 


2.5 ± 0.5 


35% 


3% 


8% 


0.640* 


0.005 


20,188 


3.733 


6.94 


7.69 


2.7 ± 0.3 


37% 


4% 


12% 


Table 1 



Sample density (p), density fluctuations (a), number of spheres in the central region 
(Ng), average number of neighbours (N c ), topological density (Tp), percent of local 
configurations with a given local order: disordered (dis), fee - like and hep -like. The 
symbols '*' in the first column indicate measures on the smaller samples (~ 35,000 
beds). The three quantities and Nc 3 \ correspond respectively to the 

average number of neighbours within the three radial distances: (1) d, (2) d+v/2 
and (3) d+v, where v is the voxel-size. 

Central Region: All the analyses reported hereafter have been performed 
over a central region (G) at 4 sphere-diameters away from the sample bound- 
aries. Note that spheres outside G are considered when computing the neigh- 
bouring environment of spheres in G. The two large samples have about 
Nq ~ 80, 000 spheres in G, whereas the four smaller have about Nq ~ 20, 000. 
In Table 1 the number of spheres in this region (Nq) is reported for each sam- 
ple. 



3 Study of the Packing Configurations and Discussion 

Densities: We calculate the local- densities and the sample-densities. The 
local-densities are the fractions between the sphere-volumes and the volumes 
of the Voronoi cells constructed around the centre of each sphere in the sample 
(the Voronoi cell is the portion of space closest to a given centre in respect 
of any other centre). The sample-densities are fractions between the sum over 
the volumes of the spheres in G and the sum over the volumes of the Voronoi 
cells associated with these spheres. We observe bell-shaped distributions with 
average densities in the range 0.586 < p < 0.640 and standard deviations a 
within 1.5 % (see Table 1). However, the local densities are not homogeneously 
distributed in the sample. Typically, the densities are relatively smaller than 
the average in a region close to the cylinders central axis; the density increases 
going out from the centre then it saturates to rather homogeneous values up 
to a distance of a few (2-3) sphere diameters from the boundary. Rather inho- 
mogeneous densities are also observed in the vertical direction, but in this case 
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Fig. 2. Packing Efficiency e v.s. radial distance (in bead diameters), (left) The 
thick lines refer to the average value within a central region; the thin lines above 
and below represent upper and lower limits for each sample; the dotted line is the fee 
configuration. The colours (online only) are associated with the different samples, 
(right) The same plot of only the average values in log-log scale; the two straight 
lines are references representing respectively the power laws: ear 2 and e oc r 4 . 

we find different behaviours depending on the sample-preparation. A detailed 
analysis of these behaviours will be the subject of a future publication. 

Local Environment: The identification of touching spheres is, in general, 
an ill-defined problem from a experimental point of view. Indeed, the result is 
affected by the precision of the calculated sphere centres and their diameters. 
At the present we are developing a new technique to overcome this problem. 
With the present data, we assume that the spheres in contact are located at 
a radial distance between the diameter d and d + v, where v is the voxel-size. 
Table 1 reports the values of the average number of contacting neighbours 
(N c ) computed in G at the three different radial distances: d, d + v/2 and 
d + v. We observe values between N c ~ 3 and 8, and an increasing trend with 
the packing density. 

Packing Efficiency: We compute the number of spheres placed within a 
certain radial distance from a given sphere. This is a measure of how effi- 
ciently local dense agglomerate of spheres are formed and therefore we called 
it packing efficiency e. It is well known that no more than 12 spheres can be 
found in contact with one sphere (the 'Kissing number' [1]), but the upper 
limit for the number of spheres within a given radial distance is, in general, 
unknown. Here, we empirically investigate the local environment surrounding 
each sphere starting from the spheres in contact up to the ones at a radial 
distance of two bead diameters. Fig. 2 shows the average, the maximum and 
the minimum numbers of neighbours within a given radial distance from any 
sphere in G. Interestingly, we observe that, in the regions around r ~ 1.4 
and 1.7, disordered packings show better average packing efficiency than fee. 
Moreover, we observe that in all the samples the upper limits are above or 
equal to the fee packing-efficiency in a large range of the radial distances. Let 
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Fig. 3. (left) Radial distribution function g(r) v.s. radial distance r (in bead diam- 
eters) for all the samples examined. The insert shows the two peaks at r = v3 and 
r = 2 in detail, (right) In the region between 1 < r < 1.4, the radial distribution 
function decreases as a power law: g(r) oc (r — ro) _a (linear trend). The insert shows 
the best-fit estimations for the exponent a, demonstrating an increasing trend with 
the density. 

us stress that, these surprisingly high efficiencies in disordered packings could 
be the inner driving cause for the occurence of such disordered assemblies. 

Radial Distribution Function: The radial distribution function g(r) (shown 
in Fig. 3) is calculated as the average number of sphere centres, within a radial 
distance r — A/2 and r + A/2, divided by cr 2 . The constant c is fixed by impos- 
ing that asymptotically g(r) — > 1 for r — > oo. (We have verified that different 
choices of A within a broad range of 10~ 4 to 10 -2 (in bead diameters units), 
lead to almost indistinguishable results.) The detail of the close-neighbour 
region (Fig. 3) shows that the two peaks at r = y/3 and r = 2 (in bead diam- 
eters) both increase in height with the packing density. This is an indication 
that there is an increasing number of configurations in the contact network [9] 
with edge-sharing in-plane triangles (peak at a/3) and with three (or more) 
aligned spheres (peak at 2). This might indicate an increasing organisation in 
the packing structure but, on the other hand, no clear signs of crystallisation 
were detected (see below). 

For all the samples investigated, we found that in the range r < r < 1.4, 
the radial distributions follow a power law behaviour: g(r) oc (r — r )~ a , with 
exponents a between 0.3 and 0.45 and the singularity at r = 1.03 (Fig.3). 
A similar behaviour, but with a = 0.5 and ro = 1, was reported in [10] for 
molecular dynamics simulations. In Fig.3 it is also highlighted the growing 
trend of the exponent a with the density p. To our knowledge, this is the 
first time that such behaviours (still theoretically unexplained) are observed 
in experiments. 

Shell Analysis: The sphere neighbourhood environment has been also stud- 
ied topologically by performing a shell analysis on the contact network and 
calculating the topological density (T D ) [11,12,13,14,15]. We observe that the 
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number of spheres at a given topological distance (j [16]) from a central one 
follows quite accurately the quadratic law: Kj = STpj 2 + c±j + c . The best 
fit values of Tp's are reported in Table 1. (This analysis was performed by 
defining in contact spheres at radial distance within r < d + v/2; the depen- 
dence of T D on the choice of the contact-neighbour criteria will be discussed 
in a forthcoming paper). We note that the topological densities increase with 
the sample-densities and are consistently smaller than 10/3 = 3.33..., which is 
the theoretical lower limit for Barlow packings [4] . A view of the topological 
shell structure is shown in Fig.l. 

Local Order: We estimate the local packing orientation by associating a 
set of spherical harmonics to the vectors r between a sphere and the set of 
neighbouring spheres (|f| < 1.2 - in beads diameters): Yi jm (0(f), 4>{r)) (with 
6(f) and <f>(f) the polar and azimuthal angles associated with r). Following 
[17], a set of rotationally invariant quantities are calculated: 



4-7T ' 

G' = \ 2ZTI £ lEnmW*), 

\ ' m=—l i£G 



We observe that the Qi's (with I — 4, 6) calculated over many sub- volumes of 
the samples tend to vanish as the sub-volumes become larger. This indicates 
that these packings are macroscopically like isotropic liquids [17] and across 
the entire system there is no prevailing orientational order. On the other hand, 
even if globally there are no signs of order or organisation, this method is useful 
to investigate the orientational properties of the local configurations. In par- 
ticular, for each sphere in the central region we measured the coefficients Q± 
and Qq calculated from the angles between a sphere and all its neighbour up 
to a radial distance r = 1.2. We first search for signatures of local crystalline 
orientational order measuring the fraction of local configurations with (Q±, Q 6 ) 
in a region close to the fee (0.1909, 0.5745) (/cc) and hep (0.0972, 0.4848) (/lcp) 
order. We observe that a relatively small fraction of configurations have 
and Qq within a range ±0.05 of their ideal fee and hep values (see £ /cc' an d 
'hep' in Table 1). We also verify that other special configurations like Icosahe- 
dral, bec and Simple Cubic, have significantly lower occurrences. On the other 
hand, we observe that large percentages (between 27% and 37%) of the lo- 
cal configurations share a common signature in their orientational order with 
0.15 < Q4 < 0.25 and 0.4 < Q 6 < 0.5 (see 'cfc' in Table 1). To understand 
whether this is an indication of a 'typical' local disordered packing around 
(0.20, 0.45)^) is the matter of present investigations. 

Structural Arrest: We study a quantity which is relevant for system dynam- 
ics: the escape probability which is the probability that a sphere can moves 
outward from a given local configuration without re- adjusting the positions 
of its first neighbours [18]. This quantity is calculated by constructing circles 
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through the centres of the three spheres corresponding to the three faces in- 
cident at each vertex in the Vorono'i polyhedron. If one of these circles has a 
radius larger than d, it implies that the central sphere can pass through that 
neighbouring configuration and move outward from its local position without 
displacing the first-neighbours. In other words, the neighbouring 'cage' is open 
if at least one radius is lager than d; vice versa the cage is closed when all radii 
are smaller than d. The escape probability is defined as the fraction of open 
cages. We find that all the samples with p > 0.6 have zero escape probability 
in G. Whereas, the two samples with p = 0.586 and p = 0.595, have very 
small fractions of open cages (0.1% and 0.6%, respectively). 

Following a theoretical approach recently developed by one of the authors [18], 
the packing realization can be considered as an inherent structure and the 
(thermo) dynamical approach toward this configuration can be reconstructed 
by virtually decreasing the sphere-diameters reducing the effective density. 
(Here, the main underlying assumption is that the system dynamics before the 
structural arrest develops through configurations around the final inherent- 
structure realization.) We find that, for all the samples with p > 0.6, the 
escape probabilities go to zero at effective densities between 0.61 and 0.63. 
This strongly suggests that around p ~ 0.62 ± 0.01 an important phase in 
the system dynamics reaches an end: above this density, local readjustments 
involving only the displacement of a single sphere are forbidden and the sys- 
tem compaction can proceed only by involving the collective and correlated 
readjustment of larger set of spheres. 
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